// Graph coloring 
#include <ctime> 
#include <stdio.h>  
#include <stdlib.h>  
#include <time.h>  
#include <iostream> 
#include <math.h>  
#include <set> 
#include <assert.h>
#include <fstream>
#include <string>
#include "graphColoring.h" 

using namespace std;  

//----------------------- Utility functions -----------------------//

int findMultiple(int multipleOf, int x) {
	int base = multipleOf;
	int powerIndex = 0;
	int num = 0;
	
	while (num <= x) {
		powerIndex++;
		num = base * powerIndex;	
	}
	
	cout << "Closest multiple of " << base << ": " << num << endl;
	return num;
}

//----------------------- Graph initializations -----------------------//

// Adapted from Grossnet 2011
// reads a matrix market format (.mtx) into an adjacency list
void readGraphAdj(int *&adjacencyList, const char *filename, int _gridSize, int _blockSize, int &graphSizeRead, int &graphSize, long &edgeSize, int &maxDegree) {		
	char comments[512];
	int graphSizeX, graphSizeY, numUsefulEdges, from, to, weightedGraph;
	int *edgesCount, *degreeList;
	double weight;
	
	numUsefulEdges = maxDegree = 0;
		
	ifstream graphFile(filename);
	
	if (graphFile.is_open()) {
		//
		// Pass 1: get the degree of each node
		//
		while (graphFile.peek() == '%') {
			graphFile.getline(comments,512);
		}
		
		graphFile >> graphSizeX >> graphSizeY >> edgeSize;
		cout << "Rows: " << graphSizeX << "  ,  Columns: " <<  graphSizeY << "   - Number of edges: " << edgeSize << endl;
		
		if (graphSizeX != graphSizeY) {
			cout << "Non Symmetric graph!" << endl;
			exit(1);
		}
		else {
			//cout << "Is it a weighted graph(1:yes  -  0:no): ";
			//cin >> weightedGraph;
			
			string S("auto.mtx");
			if(S.compare(filename) == 0) {
				weightedGraph = 0;
			}
			else { 
				weightedGraph = 1;
			}
			
			graphSizeRead = graphSizeX;
			graphSize = findMultiple(_gridSize*_blockSize, graphSizeRead);
			
			degreeList = new int[graphSize];
			memset(degreeList, 0, graphSize*sizeof(int));
			
			
			for (int i=0; i<edgeSize; i++) {
				if (weightedGraph == 1) {
					graphFile >> from >> to >> weight;	
				//	cout << from << " , " << to << " : " << weight << endl;
				}
				else {
					graphFile >> from >> to;	
				//	cout << from << " , " << to << endl;
				}
				
				
				if (from != to) { 
					degreeList[(from-1)]+= 1;
					degreeList[(to-1)]+= 1;
					
					numUsefulEdges++;
				}
			}
		}	
		
		//
		// Get max degree and reset
		//
		for (int i=0; i<graphSize; i++) {
			//cout << i << ": " << degreeList[i] << endl;  
			if (maxDegree < degreeList[i]){
				maxDegree = degreeList[i];
//				cout << i << "  Max degree: " << maxDegree << endl;  
			}
		}
	
		
		edgesCount = new int[graphSize];
		memset(edgesCount, 0, graphSize*sizeof(int));
		
		adjacencyList = new int[graphSize*maxDegree];
		memset(adjacencyList, -1, graphSize*maxDegree*sizeof(int));
		
		
		graphFile.seekg(0, ios::beg);
		
		
		
		//
		// Pass 2: Fill in the adjacency List
		//
		while (graphFile.peek() == '%') {
			graphFile.getline(comments,512);
		}
		
		
		graphFile >> graphSizeX >> graphSizeY >> edgeSize;
		//cout << graphSizeX << " " <<  graphSizeY << " - " << numEdges << endl;
		
		if (graphSizeX != graphSizeY) {
			cout << "Non Symmetric graph!" << endl;
			exit(1);
		}
		else {
			for (int i=0; i<edgeSize; i++) {
				
				if (weightedGraph == 1) {
					graphFile >> from >> to >> weight;	
					//cout << from << " , " << to << " : " << weight << endl;
				}
				else {
					graphFile >> from >> to;	
					//cout << from << " , " << to << endl;
				}
				
				
				
				if (from != to) {
					adjacencyList[(from-1)*maxDegree + edgesCount[(from-1)]] = (to-1);
					edgesCount[(from-1)]++;			
				
					adjacencyList[(to-1)*maxDegree + edgesCount[(to-1)]] = (from-1);
					edgesCount[(to-1)]++;
				}
			}
		}
	}
	else {
		cout << "Reading " << filename << " failed!" << endl;
		exit(1);
	}
	
	graphFile.close();
	
	
	edgeSize = numUsefulEdges;
	cout << "Graph: " << graphSizeRead << " - " <<  graphSize << " - " <<  edgeSize << endl;
	cout << "Max degree: " << maxDegree << endl;
	cout << "File " << filename << " was successfully read!" << endl;
	
	delete []degreeList;
	delete []edgesCount;
}


// Adapted from Grossnet 2011
// get the degree information for a graph 
int getMaxDegree(int *adjacencyMatrix, int size, int &avgDegree) {  
	int maxDegree = 0;   
	int degree;  
	avgDegree = 0;
	
	for (int i=0; i<size; i++) {  
		degree = 0;  
		
		for (int j=0; j<size; j++) {         
			if (adjacencyMatrix[i*size + j] == 1) {
				degree++; 
			}	
		}
		
		if (degree > maxDegree)  
			maxDegree = degree;  
		
		avgDegree += degree;
	}  
	
	avgDegree = avgDegree/size;
	
	return maxDegree;  
}  


// adapted from Grossnet 2011
// get the degree of each element in the graph 
void getDegreeList(int *adjacencyList, int *degreeList, int sizeGraph, int maxDegree){
    for (int i=0; i<sizeGraph; i++) {
        int count = 0;
        
        for (int j=0; j<maxDegree; j++) {
            if (adjacencyList[i*maxDegree + j] != -1)
                count++;
            else
                break;  
        }
		
        degreeList[i] = count;
    }
}


// get the nodes that are on the boundary of the subgraphs
// adapted from Grossnet 2011
int getBoundaryList(int *adjacencyList, int *boundaryList, int graphSize, int maxDegree, int _gridSize, int _blockSize, int *startPartitionList, int *endPartitionList) {  
	int boundaryCount = 0;
	set<int> boundarySet; 
	boundarySet.clear(); 
	
	int start, end, node;
	int partitionIndex = 0; 
	
	for (int i=0; i<graphSize; i++) {  		
		if (!(i < endPartitionList[partitionIndex])) {
			partitionIndex++;
		}
		
		start = startPartitionList[partitionIndex];
		end = endPartitionList[partitionIndex];
		
		for (int j=0; j<maxDegree; j++){ 
			node = adjacencyList[i*maxDegree + j]; 
			if (node != -1) {
				if ((node < start) || (node >= end)) {
					boundarySet.insert(i);
					break;		//	we just need one proof of that
				}
			}
		}
		
	} 
	
	boundaryCount = boundarySet.size();
	
	set<int>::iterator it = boundarySet.begin(); 
	for (int i=0; it != boundarySet.end(); it++) { 
		boundaryList[i] = *it;
		i++; 
	}  
	
	return boundaryCount;  
} 


//----------------------- Fast Fit Graph Coloring -----------------------//
// The sequential FF algorithm as benchmark
// adapted from Grossnet 2011
// GraphColor Adjacency list
int colorGraph_FF(int *adjacencyList, int *colors, int size, int maxDegree){  
	int numColors = 0;  
	int i, j;  
	
	int * degreeArray;  
	degreeArray = new int[maxDegree+1];  
	
	
	for (i=0; i<size; i++) {                 
		// initialize degree array  
		for (j=0; j<=maxDegree; j++) {
			degreeArray[j] = j+1;  
		}
		
		// check the colors  
		for (j=0; j<maxDegree; j++) {  
			if (i == j)	continue;  
			
			// check connected  
			if (adjacencyList[i*maxDegree + j] != -1)  {
				if (colors[adjacencyList[i*maxDegree + j]] != 0) {
					degreeArray[colors[adjacencyList[i*maxDegree + j]]-1] = 0;   // set connected spots to 0  
				}
			}
		}  
		
		for (j=0; j<=maxDegree; j++) {
			if (degreeArray[j] != 0){  
				colors[i] = degreeArray[j];  
				break;  
			}  
		}
		
		if (colors[i] > numColors) { 
			numColors=colors[i];  
		}
	}  
	
	delete[] degreeArray; 
	
	return numColors;  
}  



//----------------------- Conflict Solve -----------------------//
// adapted from Grossnet 2011
// Solves conflicts using Fast Fit natural visiting order
void conflictSolveFF(int *adjacencylist, int size, int *conflict, int conflictSize, int *graphColors, int maxDegree){
	int i, j, vertex, *colorList, *setColors;
	colorList = new int[maxDegree];
	setColors = new int[maxDegree];
	
	
	// assign colors up to maxDegree in setColors
	for (i=0; i<maxDegree; i++) {
		setColors[i] = i+1;
	}
	
	
	for (i=0; i<conflictSize; i++) {
        memcpy(colorList, setColors, maxDegree*sizeof(int));                    // set the colors in colorList to be same as setColors
		
        vertex = conflict[i]-1;
        
		
        for (j=0; j<maxDegree; j++) {                                            // cycle through the graph
			if ( adjacencylist[vertex*maxDegree + j] != -1 )                      		//      check if node is connected
				colorList[ graphColors[adjacencylist[vertex*maxDegree + j]]-1 ] = 0;
			else 
                break;    
        }
		
		
        for (j=0; j<maxDegree; j++) {                                       	// check the colorList array
			if (colorList[j] != 0){                                         //    at the first spot where we have a color not assigned
				graphColors[vertex] = colorList[j];                         //       we assign that color to the node and
                break;                                                      //   	 exit to the next
            }
        }
		
	}
}



//----------------------- Checking for error -----------------------//
// adapted from Grossnet 2011
// Checks if a graph has coloring conflicts
void checkCorrectColoring(int *adjacencyList, int *graphColors, int graphSize, int maxDegree) {
    int numErrors = 0;
	
    cout << endl << "==================" << endl << "Error checking for Graph" << endl;
	
    for (int i=0; i<graphSize; i++) {              
        int nodeColor = graphColors[i];
        int numErrorsOnRow = 0;
		
        for (int j=0; j<maxDegree;j++) {
			
            if (adjacencyList[i*maxDegree + j] == -1)
                break;
            else {     // there is a connection to that node
                int node = adjacencyList[i*maxDegree + j];
                if (graphColors[node] == nodeColor)
                {
                    numErrors++;
                    numErrorsOnRow++;
                }
            }
        }
    }
	
    cout << "Color errors for graph : " << numErrors << endl << "==================== " << endl ;   
}


//----------------------- The main driving function -----------------------//
int main(int argc, char *argv[]){  
	// usage instruction
	if (argc != 3){
		cout << "usage: " << endl 
		     << "cudaExe <input mtx file name> "
		     << " <output file name> "
		     << endl;

		return 1;
	}
	
	int maxDegree, numColorsSeq, numColorsParallel, boundaryCount, conflictCount, passes, graphSize, graphSizeRead;
	int _gridSize, _blockSize, randomnessValue;
	long numEdges;
	string inputFilename;
 
	// the adjacentList store the edge information. 
	int *adjacentList;
	
	conflictCount = 0;
	boundaryCount = 0;
	numColorsSeq = 0;
	numColorsParallel = 0; 
	
	
	//--------------------- Parameters for experiment ---------------------!

	// number of passes for GPU conflict solving	
	passes = 6;

	// the randomness from 0-2. set to 1 for the experiment
	randomnessValue = 2;


	//vertex visiting order: 0=FF, 1=SDO, 2=SL, 3=LF.
	vector<int> orderingVec;
	orderingVec.push_back(1);
	orderingVec.push_back(2);
	orderingVec.push_back(3);
		
	// recoloring iterations:
	vector<int> recoloring_iter;
	recoloring_iter.push_back(0);
	recoloring_iter.push_back(1);
	recoloring_iter.push_back(2);
	recoloring_iter.push_back(4);
	recoloring_iter.push_back(8);
		
	// synchronized recoloring? 0 for async, 1 for sync
        // Hard-coded here as we always used async
	int sync = 0;

        // recolor order
        // 0=reverse color, 1=non decreasing, 2= non increasing
	vector<int> recolor_order_vec;
	recolor_order_vec.push_back(0);
	recolor_order_vec.push_back(1);
	recolor_order_vec.push_back(2);
		

	
	cout << endl << "!--------------- Graph Coloring  -------------------!" << endl;

	// Grid and block size
	// number of threads in a cuda block. 
	_blockSize = 1024;
	cout << endl;
	
	// take input filename from the command line argument.
	inputFilename = argv[1];
	
	// array of subgrid size
	int subGraphSize[11] = {4,5,6,7,8,9,10,15,20,30,45}; //11 different sizes

	// hard code graph size for the benchmark graphs, for experimenting purpose. 
	int graphSize_a;
	if (inputFilename == "auto.mtx")
		graphSize_a = 448695;
	else if (inputFilename == "bmw3_2.mtx")
		graphSize_a = 227362;
	else if (inputFilename == "hood.mtx")
		graphSize_a = 220542;
	else if (inputFilename == "ldoor.mtx")
		graphSize_a = 952203;
	else if (inputFilename == "msdoor.mtx")
		graphSize_a = 415863;
	else if (inputFilename == "pwtk.mtx")
		graphSize_a=217918;


	//Create a file to write o/p results. the o/p filename is taken from the second command line argument
        ofstream out;
	string outputFilename = argv[2];
        out.open(outputFilename.c_str());
	
        if(!out){
    	   cerr<<"ERROR:result.txt could not be opened\n";
    	   exit(1);
        } 
    
        // output header line
	out<<"gridSize  ,"<<"#threads  ,"<<"ordering  ,"<<"#recoloring_iter  ,"<<"recoloring_order  ,    #colors,  normalized running time,"<<endl;

	// 1st level for loop: iterate through different subgraph size
	for(int i_gridSize = 0; i_gridSize < 11; i_gridSize++){
		// calculate the gridSize that based on the graph size and subgraph size. 
		_gridSize = graphSize_a/(_blockSize*subGraphSize[i_gridSize]) + 1;
		int numRandoms = _gridSize*_blockSize*10;
		
		int *startPartitionList = new int[_gridSize*_blockSize];	
	        int *endPartitionList = new int[_gridSize*_blockSize];		

	        memset(startPartitionList, -1, _gridSize*_blockSize*sizeof(int));
	        memset(endPartitionList, -1, _gridSize*_blockSize*sizeof(int));
		
		// read graph from file to adjacentlist
		readGraphAdj(adjacentList, inputFilename.c_str(), _gridSize, _blockSize, graphSizeRead, graphSize, numEdges, maxDegree);

		cout << "Allocation successful!" << endl;
	

		// initialize graphcolors, boundaryList, and degreeList array
		int *graphColors = new int[graphSize];          
		int *boundaryList = new int[graphSize]; 
		int *degreeList = new int[graphSize];
		
		memset(graphColors, 0, graphSize*sizeof(int)); 
		memset(boundaryList, 0, graphSize*sizeof(int)); 
		memset(degreeList, 0, graphSize*sizeof(int)); 
		
		
		// Get degree List
		cout << "get degree list" << endl;
		getDegreeList(adjacentList, degreeList, graphSize, maxDegree);
		

		// initialize randomList
		int *randomList = new int[numRandoms];
		for (int i=0; i<numRandoms; i++)					// stores random numbers in the range of 0 to 2
			randomList[i] = rand()%(randomnessValue+1);
		
		
		cout << "allocating exact partition size " << endl;	
		int partitionSizes = graphSize / (_gridSize*_blockSize);
	
	
		for (int i=0; i<_gridSize*_blockSize; i++){
			startPartitionList[i] = i*partitionSizes; 
			endPartitionList[i] = (i+1)*partitionSizes;
		}
		
		
		//---------------------Obtain Boundary List ---------------------!
		cudaEvent_t start_b, stop_b;
		// timer for obtaining boundary
		float elapsedTimeBoundary;
		
		cudaEventCreate(&start_b); 
		cudaEventCreate(&stop_b); 
		cudaEventRecord(start_b, 0); 
		
		cout << "obtain boundary list" << endl;
		// get boundaryCount and get boundary list
		boundaryCount = getBoundaryList(adjacentList, boundaryList, graphSize, maxDegree, _gridSize, 
										_blockSize, startPartitionList, endPartitionList);							cout << "boundary nodes count: " << boundaryCount << endl;
		
		cudaEventRecord(stop_b, 0); 
		cudaEventSynchronize(stop_b); 
		cudaEventElapsedTime(&elapsedTimeBoundary, start_b, stop_b); 
		cout << "Time to getBoundaryList :"<< elapsedTimeBoundary << " ms" << endl;
	

		
		//--------------------- Sequential Graph Coloring ---------------------!
		cudaEvent_t start, stop, stop_1, stop_4;         
		float elapsedTimeCPU, elapsedTimeGPU, elapsedTimeGPU_1, elapsedTimeGPU_4; 
		
		cudaEventCreate(&start); 
		cudaEventCreate(&stop); 
		cudaEventRecord(start, 0);  
		
		
		cout << "coloring sequentially" << endl;
		
		// color the graph sequentially using natural vertex visiting order. 
		colorGraph_FF(adjacentList, graphColors, graphSize, maxDegree);  
		

		cudaEventRecord(stop, 0); 
		cudaEventSynchronize(stop); 
		cudaEventElapsedTime(&elapsedTimeCPU, start, stop); 
		
		
		// Get colors
		numColorsSeq = 0;
		for (int i=0; i<graphSize; i++){
			if ( numColorsSeq < graphColors[i] )
				numColorsSeq = graphColors[i];
		}
		
		
		
	
		//--------------------- Checking for color conflict ---------------------!
		
		cout << endl << "Sequential Conflict check: ";
		checkCorrectColoring(adjacentList, graphColors, graphSize, maxDegree);
		
		cout << endl;  
		cout << "Sequential colors: " << numColorsSeq << endl;
	
		
	
	
	
	
	//--------------------- Parallel Graph Coloring ---------------------!	


		// 2nd level for loop: iterate through different vertex visiting ordering. 
		for(unsigned int i_ordering = 0; i_ordering < orderingVec.size(); i_ordering++){
			int ordering = orderingVec[i_ordering];
			
			// 3rd level for loop: iterate through recoloring iterations. 
			for(unsigned int i_recoloring_iter = 0; i_recoloring_iter < recoloring_iter.size(); i_recoloring_iter++){
				int recoloring = recoloring_iter[i_recoloring_iter];
				
				// if SDO order, do not recolor
				if (ordering == 1 && recoloring > 0)
					break;
				
				// 4th level for loop: iterate through recoloring order
				for(unsigned int i_recolor_order = 0; i_recolor_order < recolor_order_vec.size(); i_recolor_order++){
					int recolor_order = recolor_order_vec[i_recolor_order];
					// if recoloring iteration = 0, no need to iterate through recolor order.					
					if(recoloring == 0 && recolor_order > 0){
						break;
					}
					
					// write experiment parameters to file
					out<<_gridSize<<","<<ordering<<","<<recoloring<<","<<recolor_order<<",";

					// conflict array 
					int *conflict = new int[boundaryCount];
					
					// initialize conflict array to 0                     
					memset(conflict, 0, boundaryCount*sizeof(int));                   	    
					
					// reset colors to 0 
					memset(graphColors, 0, graphSize*sizeof(int));                         
	
	
					//--------------- Parallel Graph coloring in GPU + Conflict Detection in CPU -------------------
					cudaEventCreate(&start); 
					cudaEventCreate(&stop); 
					cudaEventCreate(&stop_1); 
					cudaEventCreate(&stop_4); 
					cudaEventRecord(start, 0); 
	
	
					int *conflictTmp = new int[boundaryCount];
					memset(conflictTmp, 0, boundaryCount*sizeof(int));  
					
					
					// here it goes the method to be run in cuda
					cudaGraphColoring(adjacentList, boundaryList, graphColors, degreeList, conflictTmp, boundaryCount, 
									  maxDegree, graphSize, passes, _gridSize*_blockSize, _gridSize, _blockSize,
									  startPartitionList, endPartitionList, randomList, numRandoms, ordering, recoloring, sync, recolor_order);
	
	
					cudaEventRecord(stop_1, 0); 
					cudaEventSynchronize(stop_1); 
		
	
					// count number of parallel colors
					int interColorsParallel = 0;
					for (int i=0; i<graphSize; i++) {
						if ( interColorsParallel < graphColors[i] ) {
							interColorsParallel = graphColors[i];
						}
					}
	
	
					//-------- Conflict Count
					conflictCount = 0;
					for (int i=0; i<boundaryCount; i++) {
						int node = conflictTmp[i];
		
						if (node != -1) {
							//cout << "conflict " << conflictCount <<  " at " << node << endl;
							conflict[conflictCount] = node;
							conflictCount++;
						}
					}
  					delete[] conflictTmp;
	
	
					cudaEventRecord(stop_4, 0); 
					cudaEventSynchronize(stop_4); 
	
					cout << "Boundary Count: " << boundaryCount << endl;
					cout << "Conflict count: " << conflictCount << endl;
					cout << endl;

					cout << "Colors before solving conflict: " << interColorsParallel << endl;
	
					//--------------- Step 4: solve conflicts 
					conflictSolveFF(adjacentList,  graphSize, conflict, conflictCount, graphColors, maxDegree); 
	
					cudaEventRecord(stop, 0); 
					cudaEventSynchronize(stop); 
	
	
					cudaEventElapsedTime(&elapsedTimeGPU, start, stop); 
					cudaEventElapsedTime(&elapsedTimeGPU_1, start, stop_1); 
					cudaEventElapsedTime(&elapsedTimeGPU_4, start, stop_4); 
	
					cout << "conflict solved" << endl;
	
					// calculate the total number of colors using GPU after conflict solving
					numColorsParallel = 0;
					for (int i=0; i<graphSize; i++){
						if ( numColorsParallel < graphColors[i] )
							numColorsParallel = graphColors[i];
					}
						
	
					//--------------------- validating colors ---------------------!
	
					cout << endl <<  "Parallel Conflict check:";	
	
					checkCorrectColoring(adjacentList, graphColors, graphSize, maxDegree);

					//--------------------- display information ---------------------!	

                	                cout << "ordering: " << ordering << "    recoloring: " << recoloring << "     recolor order: " << recolor_order << endl; 
					cout << "Grid Size: " << _gridSize << "    Block Size: " << _blockSize << "     Total number of threads: " << _gridSize*_blockSize << endl;
					cout << "Graph average subsize: " << graphSize/(_gridSize*_blockSize) << endl;
	
					cout << "GPU Passes done: " << passes << endl;
					cout << "CPU time (First Fit): " << elapsedTimeCPU << " ms    -  GPU Time (SDO Solver): " << elapsedTimeGPU << " ms" << endl;

					cout << "Sequential Colors: " << numColorsSeq << "      -       Parallel Colors: " << numColorsParallel << endl;     
					cout << "GPU normalized running time (including boundary): " << elapsedTimeCPU/(elapsedTimeBoundary + elapsedTimeGPU) << " x" << endl;
					
					//Write the result on o/p file
					out<<numColorsParallel<<","<<elapsedTimeCPU/(elapsedTimeBoundary + elapsedTimeGPU)<<endl;

					cout << "!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~!" << endl << endl;
					
					delete []conflict;
					
     			        } // end of for (recoloring ordering)
   			} // end of for (recoloring iterations)
  		} // end of for (vertex visiting order)

		// release memory
		delete []graphColors; 
		delete []startPartitionList;	
		delete []endPartitionList;
		delete []boundaryList;
		delete []degreeList;
		delete []randomList;
		delete []adjacentList;
 
	} // end of for (grid size)
		
 	out.close();	// close output file
	
	return 0;  
}  

//End of File
